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Abstract 

In the mean-field approximation, a trapped Bose-Einstein condensate at zero 
temperature is described by the Gross-Pitaevskii equation for the condensate, or 
equivalently, by the hydrodynamic equations for the number density and the cur- 
rent density. These equations receive corrections from quantum field fluctuations 
around the mean field. We calculate the semiclassical corrections to these equations 
for a general time-independent state of the condensate, extending previous work to 
include vortex states as well as the ground state. In the Thomas-Fermi limit, the 
semiclassical corrections can be taken into account by adding a local correction term 
to the Gross-Pitaevskii equation. At second order in the Thomas-Fermi expansion, 
the semiclassical corrections can be taken into account by adding local correction 
terms to the hydrodynamic equations. 
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1 Introduction 



The achievement of Bose-Einstein condensation in trapped atomic vapors ^ has rekindled 
interest in the theory of nonhomogeneous interacting Bose gases. This problem can be con- 
veniently formulated as a problem in quantum field theory, but the quantum field equations 
are extremely difficult to solve in general. Fortunately, many of the basic properties of the 
condensates in existing experiments, such as the density profile of the ground state and 
the frequencies for small amplitude collective oscillations of the condensate, can be de- 
scribed reasonably well using the mean-field approximation. This approximation reduces 
the problem to solving the Gross-Pitaevskii equation, a classical equation for the mean 
field. However, fluctuations of the quantum field around the mean field provide corrections 
to mean-field predictions that grow as the square root of the number density of the atoms. 
These corrections will become increasingly important as higher condensate densities are 
achieved and as the precision of experimental measurements improves. Furthermore, there 
are some observables, such as the damping rates of collective oscillations, that vanish in the 
mean-field approximation and are therefore sensitive to the effects of quantum field fluc- 
tuations. It is therefore important to understand the effects of quantum field fluctuations, 
and to be able to calculate them accurately. 

The problem of the nonhomogeneous interacting Bose gas at low temperature was 
studied by Fetter in 1972 using the Bogoliubov approximation. In this approximation, 
the quantum field equation is linearized around the solution to the Gross-Pitaevskii 
equation. More elaborate approximations that provide a better extrapolation to higher 
temperatures, such as the Hartree-Fock-Bogoliubov (HFB) approximation, were developed 
in the early 1980's p, |[. They require the self-consistent solution of a coupled system of 
equations consisting of a partial differential equation for (p and a linear equation for the 
quantum fluctuation field ip, both of which involve expectation values of operators quadratic 
in -0. A critical analysis of the HFB approximation and other simpler approximations has 
recently been given by Griffin 

One of the attractive features of the HFB approximation is that it is a conserving 
approximation that is guaranteed to respect the conservation laws that follow from the 
symmetries of the field theory. One of the problems with the HFB approximation is that it 
is not a gapless approximation. In the case of a homogeneous gas, the HFB approximation 
gives a gap in the spectrum in violation of Goldstone's theorem. More generally, it fails 
to respect some of the consequences of the spontaneously broken symmetry of the theory. 
The classification of approximation methods for a Bose gas according to whether they are 
conserving or gapless approximations was first made by Hohenberg and Martin 0. A 
controlled approximation that corresponds to the truncation of a systematic expansion in 
a small parameter is guaranteed to be gapless. However such approximations are typically 
not conserving, because the conservation laws may be satisfied only up to corrections that 
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are higher order in the expansion parameter. 

Another problem with the HFB approximation is that solving the HFB equations is a 
complicated computational problem. It requires solving the partial differential equation for 
(j), solving the eigenvalue equations for all the normal modes of tp, calculating expectation 
values of operators quadratic in ip, and iterating this sequence of calculations until they are 
self-consistent. The resulting solutions for and for the normal modes of contain a great 
deal of information. Unfortunately, because the HFB approximation is not gapless, some 
of that information is qualitatively wrong. Furthermore, much of that information is not 
easily accessible in experiments. For example, the total contribution to the number density 
from all the normal modes is easier to measure than the contribution from individual normal 
modes. 

It would be worthwhile to develop an improved treatment of the nonhomogeneous 
interacting Bose gas that avoids the inconsistencies of the HFB approximation and is 
also computationally simpler. It was shown in Ref. that the effects of quantum field 
fluctuations on the ground state can be taken into account by adding local correction terms 
to the partial differential equations of the mean-field approximation. The consistency of 
the approach was guaranteed by using only controlled approximations. The enormous 
simplification compared to the HFB approximation comes at the expense of any direct 
information on the individual normal modes of the quantum fluctuation field. In this 
paper, we extend the results of Ref. [0 to arbitrary time-independent states of the Bose- 
Einstein condensate at zero temperature, including vortex states. We also streamline the 
derivations in Ref. by using dimensional regularization to control infrared and ultraviolet 
divergences. 

The controlled approximations that guarantee the consistency of our approach are 
based on truncations of systematic expansions in two small quantities. The first of these 
quantities is [64p(r)a^/'7r]^/^, where a is the S-wave scattering length and p is the number 
density. This quantity is a local measure of the magnitude of the effects of quantum field 
fluctuations. We restrict our analysis to the semiclassical approximation, in which the 
expansion is truncated at first order in y/ pa'-^. The second small quantity is ^/-R, where 
^ = [167rap(r)]~^/^ is the local coherence length and R is the length scale for significant 
variations in p(r). The expansion in ^/i? defines the Thomas- Fermi approximation. The 
effects of quantum field fluctuations are dominated by modes with wavelengths of order 
^. The quantity C,/R is therefore a measure of how local the effects of the quantum field 
fluctuations are. 

In present experiments with Bose-Einstein condensates, the largest values of the semi- 
classical expansion parameter have been achieved for condensates of Na^^. The scattering 
length of Na^^ is a = 2.75 nm. In the experiment of Ref. [§], the peak density that was 
achieved was approximately p = 3 x 10^^ atoms/cm^. The peak value of the semiclassical 
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expansion parameter is therefore (QApa? /txY^"^ = 0.04. This is small but still large enough 
that semiclassical corrections should be observable in quantities that can be measured with 
high precision, such as the collective excitation frequencies of the condensate. The value of 
the coherence length ^ at the peak density is [167rap(r)]~^/^ = 5 x 10~^ cm. The transverse 
radius R of the elongated condensate was roughly 3 x 10^'^ cm. Thus the Thomas-Fermi 
expansion parameter ^/R a,t the peak density is roughly 0.02. The sizes of the expansion 
parameters can be controlled by adjusting the shape of the trapping potential and the 
number of atoms in the trap. The number density that can be achieved is limited however 
by the loss of atoms from the trap due to 3-body collisions. The size of the expansion 
parameters can also be controlled by changing the scattering length of the atoms. For 
example, the use of Feshbach resonances to change the scattering length of Na^^ atoms has 
already been demonstrated 0. 

In the mean-field approximation, the condensate for a time-independent state of 
a Bose- Einstein condensate in an external potential V^(r) satisfies the Gross-Pitaevskii 
equation: 

= 7^V2 + /i-yU 1 

\2m j m 

where a is the S-wave scattering length and /i is the chemical potential. This equation 
can be written equivalently in the hydrodynamic form that consists of a pair of coupled 
partial differential equations for the number density p = |0p and the current density 
j = _^(0*V0-(/)V0*): 

^i = V + p + -2pVV+ Vpf + -^f, 2 

m 8mp^ L J 2p^ 

= V-j. (3) 



We calculate the corrections to the mean-field equations from quantum field fluctua- 
tions using the semiclassical approximation, which includes all terms through first order 
in y/pa^. We use the Thomas-Fermi expansion to express the corrections in terms of local 
quantities. In the Thomas- Fermi limit, we find that the quantum field fluctuations can be 
taken into account by adding a local correction term to the Gross-Pitaevskii equation (|I|): 



2m 



V^ + p-V 



Anh a, 



m 




(4) 



At second order in the Thomas- Fermi expansion, the correction term is no longer local when 
expressed in terms of 0. However there are also nonlocal terms in the relation between 
p and (f). We find that the nonlocal terms cancel to give a local correction term to the 
hydrodynamic form of the mean-field equations. The current density j still satisfies the 
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continuity equation (H), while (|D is replaced by 



P = V + p + -2pV^p+{Vpf + — j2 

m Qmp'^ L -I ip'^ 

This equation was derived previously in Ref. ||^ for the special case of the ground state for 
which j vanishes. In this paper, we extend that derivation to a general time-independent 
state with a nonvanishing current. We find a remarkable cancellation of the j-dependent 
terms in the semiclassical correction, so that the dependence on the current j enters only 
through the classical terms. 

We begin in Section ^ by writing down the general equations that determine the con- 
densate, the number density, and the current density for a time-independent state of a 
Bose-Einstein condensate in an external potential. In Section |^, we set up a perturbative 
framework for calculating the local effects of quantum field fluctuations when they are 
dominated by wavelengths comparable to the coherence length. The semiclassical approxi- 
mation is introduced in Section ^ The equation for the condensate in this approximation is 
expressed in terms of p, j, and expectation values involving the quantum fluctuation field. 
The Thomas-Fermi expansion is used in Section ^ to calculate the expectation values to 
second order in the gradient expansion. Inserting these expectation values into the semi- 
classical equation, we obtain our final result (^). In Section ^, we deduce the corresponding 
equation for the condensate and show that it is logarithmicly sensitive to wavelengths 
much larger than the coherence length. Finally, in Section |^ we summarize our results 
and discuss some possible generahzations. Calculational details related to the evaluation 
of Feynman diagrams are collected in two appendices. In Appendix 0, we give expressions 
for momentum integrals with dimensional regularization used for both an infrared and an 
ultraviolet cutoff. In Appendix |B|, we illustrate the use of Feynman diagrams to calculate 
the coefficients in the Thomas- Fermi expansions of the expectation values of local quantum 
field operators. The coefficients that appear in the equations for the mean field and for 
the hydrodynamic variables are given in terms of dimensionally regularized integrals in 
Appendix |C[ 



2 Static Bose-Einstein condensate 



The problem of the Bose-Einstein condensation of a large number of identical atoms 
in a trapping potential V{r) can be formulated in terms of a quantum field theory with a 
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single complex field ipiv^t). The quantum field satisfies 



+ + /i - V{r) ip 



0. 



(6) 



For simplicity of notation, we have set h = 2m = 1 in (^. Dimensional analysis can be 
used to reinsert the appropriate factors of h and 2m at the end of the calculation. The 
parameter g in (^) is proportional to the S-wave scattering length a of the atom: 

g = = IGvra. (7) 

m 

The chemical potential /j, in is to be adjusted so that the average number of atoms in 
the ground state |0) is a specified number A^: 



N. 



(8) 



We now consider an arbitrary time-independent state \n) of the Bose-Einstein conden- 
sate at zero temperature. If the potential traps the atoms in a simply-connected region of 
space, the time-independent states are the ground state and vortex states. In the ground 
state, the condensate 0(r) has a constant phase and the current density j(r) vanishes. In 
a vortex state, j is nonzero and the phase of changes by a multiple of 2tt along a curve 
that circumscribes the core of the vortex. If the atoms are trapped in a region that is not 
simply-connected, there can be more complicated states |n) with time-independent number 
density and current. We denote the number density and the current density in the state 
|n) by p and j, respectively: 



p = {n\^''ijj\n), 



(9) 
(10) 



We also denote the condensate, which is the expectation value of ip in that state, by 0: 

= {n\ip\n). (11) 

The quantum fluctuation field ■0, which is defined by 

^(r,t) = 0(r) + Vi(r,t), (12) 

has a vanishing expectation value: {n\ilj\n) = 0. The quantum field equation (^) can be 
expressed as a classical field equation for coupled to a quantum field equation for -0: 

9_ 
2 








iijj + [V^ + - 1/ -ijj 



20 + 2(^t^)0 + (V^^)0* + {^^tPi^) 



(13) 



2^ + 02^t + 



+0* [ij^j - ii^ij)^ + (vV^ - (^VV") 



(14) 
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The equation ([T3|) for is just the expectation value of (|^). The equation (|^) for ip is 
obtained by subtracting ( ]TB| ) from (^). We have suppressed the state |n) in the expectation 
values that appear in ( |I3| ) and (|14|). Through the remainder of this paper, it will be 
understood that all expectation values are for the time-independent state \n). 



The system of equations consisting of (|T^) and (|T4D involves coupled nonlinear equa- 
tions for (f) and and is extremely difficult to solve in general. We must resort to some 
approximations to make this system of equations more tractable. The consistency of 
an approximation will be guaranteed if it is controlled by a small expansion parameter. 
From the example of the homogeneous Bose gas, it is clear that the fractional size of 
the effects of quantum field fiuctuations is measured locally by the dimensionless quantity 
[64p(r)a^/7r]^/^. If the peak value of this quantity is less than one, one can use y/pa^ as 
an expansion parameter. If we wish to describe the interactions of the atoms in terms of 



the S-wave scattering length a only, then we are limited to first order in 



At second 



order in 



there are ultraviolet divergences whose renormalization requires additional 



input parameters, including a pointlike contribution to the 3^3 scattering amplitude 
rri [1^ . We will refer to the approximation in which we truncate at first order in i/po^ as 



the semiclassical approximation. 

The semiclassical approximation is defined not by neglecting any particular terms in 

^pa^ at first 



the equations (|T3[) and (0), but rather by the truncation of the expansion in 
order. This approximation does allow those equations to be simplified, because not all the 
terms are necessary to achieve first order accuracy, but the simplifications depend on the 
quantity being calculated. Our goal will be to achieve first order accuracy in quantities 
like 0, p, and j that can be defined as expectation values of local operators. In terms of 
Feynman diagrams, first order accuracy requires keeping all contributions from one-loop 
diagrams, while diagrams with 2 or more loops can be neglected. The expectation value 
{ip'^ ipip) receives contributions only from diagrams with 2 or more loops, and it therefore 
can be omitted in (^) . We need to solve the quantum field equation (|T^ only with enough 
accuracy to guarantee first order accuracy for the remaining expectation values (ip^ip) and 
{ipip) in (HID- The terms in ([I^) that are quadratic or cubic in ip contribute to (ip'^ip) 
(ipip) only through diagrams with 2 or more loops. We therefore need keep only those 
terms in (14) that are linear in ip. Thus to the accuracy required, (|l^) and (14) reduce to 











iip + + p-V ip ~ gm ip 



||0|20 - g{i,^i,)<P - 

2„ 



(15) 
(16) 



The expectation value {'ip'^ip) in (0) is the noncondensate density, which is the correction to 
the number density from quantum field fiuctuations. The expectation value (ipip) is called 
the anomalous density. It should be emphasized that the equation (|16D is not designed to 
give the normal modes of ip to any specific accuracy. We need only solve the equation for 
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ijj with sufficient accuracy so that the expectation values in (|T5D can be calculated with 
an error that is second order in y/ pa'-^ relative to Since these expectation values are 
dominated by modes with wavelengths comparable to the coherence length ^, the equation 
for need only be accurate for wavelengths of order ^. 

We now compare the accuracy of the semiclassical equations (^) and ([T6|) to the tra- 
ditional approximations to the full equations (0) and (|14D . While the main purpose of the 
traditional approximations was to take into account the effects of nonzero temperature, it 
is instructive to determine their accuracy at zero temperature. The Bogoliubov approxi- 
mation is defined by the equations (^) and (|TBp, with the expectation values set to zero 
in ([T5|). Thus (j) satisfies the Gross-Pitaevskii equation. The solution to this equation gives 
0(r) with a fractional error that is first order in y/pa?. In contrast, the solution to the 
semiclassical equations (|1^) and ([16|) gives (p with a fractional error that is second order in 



In the Hartree-Fock-Bogoliubov (HFB) approximation, the {ip'^ipip) term in ( p!3D is 
neglected so that it reduces to (|T5|) and the products of field operators in (|l^) are replaced 



by terms in which pairs of operators have been contracted into expectation values: 







iip + + p-V ^ 



(17) 



Solving the system of equations ([T5| ) and ( pT]) is a complicated numerical problem. It 
requires making an initial guess for </> and for the infinitely many normal modes of ip, 
calculating the expectation values (ip^ip) and (ipip), solving (|T5|) for (p, solving (p!7| ) for 
the normal modes of ip, and iterating this sequence of calculations until they are self- 
consistent. The net effect of the additional terms in (p!7| ) that do not appear in (|T6D is to 
include a subset of corrections that are second order or higher in y/pa?. Since there are 
other corrections that are second order in \/pa? that have not been included, the use of ( [T^ ) 
in conjunction with ([l5|) does not improve the accuracy in the determination of (p. If all 
the second order corrections were uniformly small in y/pa^, it would do no harm to include 
only some of them. Unfortunately, there are individual correction terms that diverge in 
the limit ^/i? —>■ 0, but which cancel when all terms of a given order in -y/po^ are added 
together. Including some but not all of the second order terms is dangerous, because one 
may omit some of the terms that are necessary for the cancellations. Thus, by including a 
subset of the higher order corrections, the HFB approximation may actually decrease the 
accuracy. This problem can be avoided by using a controlled approximation, because any 
such cancellations will occur order by order in the expansion parameter. 

The fact that the HFB approximation is not a controlled approximation is reflected at 
nonzero temperature in the existence of an energy gap in the spectrum of the hamiltonian 
for ip, which contradicts Goldstone's theorem. In the Popov approximation, the energy gap 
is ehminated by setting the anomalous density {'ipip) to zero in (|l^) and (|l^). The resulting 
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error in is first order in \/pa^. Thus at zero temperature, the Popov approximation is no 
more accurate than the Bogohubov approximation, despite the fact that it requires much 
more complicated calculations. In the approximation used by Goldman et al. in Ref. 0], 
the {ipi^) term in (|15D was set to zero. Thus it also gives an error in cj) that is first order in 



The semiclassical equations (^) and ([l6|) can be greatly simplified by using a further 
approximation that is also controlled by a small expansion parameter. The effects of 
quantum field fluctuations are dominated by wavelengths comparable to the local coherence 
length ^ = [167rp(r)a]~^/^. For the ground state of a Bose-Einstein condensate containing 
atoms, ^ scales like N~'^^^, where N is the number of atoms in the trap [|10|- In contrast, 
the length scale R for significant variations in p(r) grows like N^^^. If is sufficiently large, 
^ is much shorter than R. This justifies an expansion in powers of ^/i? or, equivalently, 
in powers of gradients of 0. This expansion defines the Thomas-Fermi approximation. 
We will calculate the semiclassical corrections through second order in the Thomas-Fermi 
expansion. 

In the Thomas- Fermi limit, the semiclassical equation (|T3|) for the condensate becomes 
particularly simple. As we shall see, the expectation values at leading order in C,/R are 
simply functions of 0: 



1 



Inserting the expectation values (|18|) and (|19|) into (]T5|), it reduces to a partial differential 
equation for 0. Using dimensional analysis to reinsert the appropriate factors of 2m and h, 
we obtain the equation (^, which is simply the Gross-Pitaevskii equation with a correction 
term that corresponds to the local density approximation. 

To go to second order in the Thomas- Fermi expansion, we need to include terms in the 
expectation values ( ]T8| ) and (|19|) that are second order in gradients of 0. As we shall see, 
the coefficients of these terms are logarithmicly infrared divergent, indicating a sensitivity 
to length scales much greater than ^. This sensitivity can be reduced by expressing the 
equation (|15D in terms of the hydrodynamic variables p and j defined by @ and (0): 



p = |0|2 + (20) 
j = _i (0*V0 - 0V0*) - -Vi!^). (21) 

We will find that if the condensate in (|15]) is eliminated in favor of p and j, the infrared 
divergences cancel and the equation reduces to (|^). We will obtain this result in two 
different ways. In Sections § and |^, we organize the calculation so that ([TsD gives directly 
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our final equation for p and j. In Section ^ we calculate the expectation values (|T^) 
and ( [TP| ) to second order in the Thomas-Fermi expansion using an infrared cutoff. We then 
use the condensate equation ( ]T5| ) to deduce our final equation for p and j. 



3 Perturbative framework 



In this section, we set up a perturbative framework that can be used to calculate the 
leading effects of quantum field fluctuations. This framework involves introducing an arbi- 
trary momentum scale A that will later be chosen to be the inverse of the local coherence 
length ^. The perturbation series is an expansion in powers of (?A, and the semiclassical 
approximation introduced in Section ^ is the truncation of that expansion at first order in 

The quantum field theory associated with the equation is summarized by the action 

S[ilj] = J dt J d^xS^^^iP^-ij^tp) -V^^ -V^ + [p-V{r)]4j^4j - (22) 

This action is invariant under the U{1) symmetry in which the field is multiplied by a 
phase: 

^{v,t) e^Xr,t). (23) 
The action is also invariant under time-reversal symmetry: 

ij{r,t) — . ^^(r,-t). (24) 



The U{1) symmetry (pSf ) guarantees that the number density and the current density 
satisfy the continuity equation: 

P + V-j = 0. (25) 

A U{1) symmetry transformation can also be used to make the condensate 0(r) = (ipir)) 
real- valued at any specific point tq. 

The quantum field theory defined by (^) has ultraviolet divergences that must be 
removed by renormalizations of the parameters p and g. There is also an ultraviolet di- 
vergence in the number density that can be removed either by renormalization or by an 
operator-ordering prescription. In Ref . [0 , the effects of quantum field fluctuations on the 
ground state were calculated using a momentum cutoff to regularize the ultraviolet diver- 
gences. The divergences were cancelled explicitly by counterterms for p, g, and p. In this 
paper, we choose to regularize the ultraviolet divergences using dimensional regularization, 
which involves calculating integrals as analytic functions of the number of dimensions D 
and then analytically continuing them to = 3. This method has been used to calculate 
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the properties of a homogeneous Bose gas [ITl, O] and it greatly streamhnes the calcula- 



tions. One advantage of dimensional regularization is that it sets integrals that diverge 
like a power of the ultraviolet cutoff, and involve no other momentum scales, equal to zero. 
Since the counterterms for /x, g, and p are pure power ultraviolet divergences, they are 
identically zero in dimensional regularization. We have therefore omitted counterterms for 
the parameters g and p in (p2D . 



It is convenient to decompose both the condensate (|T^) and the quantum fluctuation 
field ip defined by (|12D into real and imaginary parts: 

0(r) = i;(r)+2w(r), (26) 
^.(M) = ^'■••'>^"'<^-') . (27) 

In the neighborhood of a point tq where w vanishes, we can identify ^ and r] as the quantum 
fluctuation fields associated with the number density and with the phase of the condensate, 
respectively. Inserting these expressions into the action (|2^) and expanding in powers of ^ 
and ?7, we obtain a condensate term S[v + iw\ and terms up to fourth order in ^ and rj. The 
equations (p!2[) and (^) define the Cartesian parameterization of the quantum field if) in 
terms of real quantum fields. One could equally well use an alternative parameterization, 
such as the polar parameterization: 

^(r,t) = [t;2(r) + a(r,t)]^^%*^We'"('"'*\ (28) 

where v and 9 are real classical fields, while a and a are real quantum fields. The polar 
parameterization has the advantage in perturbative calculations that it eliminates infrared 
divergences that with other parameterizations cancel only after adding Feynman diagrams 
||14||. If we use controlled approximations based on systematic expansions in small quanti- 



ties, any parameterization should lead to the same results for physical quantities. In Ref. 
, it was verified explicitly that the Cartesian and polar parameterizations give the same 
final equation for the number density in the ground state. Since the calculations were 
simpler using the Cartesian parameterization, we use only that parameterization in this 
paper. 

If the system is near the Thomas- Fermi limit, the length scale for significant changes 
in the condensate is large compared to the local coherence length [gp]^^^"^. At any specific 
point To inside the condensate, the short-wavelength modes of the quantum fields ^ and 
?7 behave locally like those of a homogeneous Bose gas with coherence length [gpivQ)]'^^"^ . 
Their dispersion relation can be approximated by the Bogoliubov dispersion relation 



e(fc) = + A2, (29) 

with ^ gpiro). A particle with the dispersion relation (129|) is described by a free field 
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theory with the action 



dt / rf^x <{ i (r]^ - ^T] 



(V0' + AT + (Vr7) 



(30) 



We can express the complete action ( p2D as the sum of the background term S[v + iw], the 
free action (BO) for the quantum fields, and an interaction term: 



S[ip] = S[v + iw] + %ee[^.^] + Si^-,t[v,w,^,r]]. 



(31) 



The last term describes the interaction of the Bogoliubov modes with each other and with 
sources that depend on v and w: 



j dtj d^x |V2(Te + St]) + ^(Xe' + Yrf + 2/?^^?) 



+ —[Zi + Qr^){e + V^) 



The sources are 



T 
S 
X 
Y 
Z 

Q 

R 



/i-l^(r)-^(t;2 + ^i;2) 



A2 + ^_y(r)_|(3t;2 + ^2^, 

^Ji-V{v)-^-{v'' + ?>w^), 
9 

V, 

2 ' 
9 

—gvw. 



2 



w + 'V w, 



(32) 

(33) 
(34) 
(35) 
(36) 
(37) 

(38) 
(39) 



The sources T and S are simply the real and imaginary parts of the Gross-Pitaevskii 
equation for the condensate v + iw. Note that the dependence of the source X on 
precisely cancels the A-dependence of the free action (|30|). Thus A is a completely arbitrary 
parameter. That arbitrariness will be exploited below to ensure that the terms in (|3^) can 
be treated as perturbations. 



The quantum field equations for ^ and r] are obtained by varying the action 
Taking the expectation value of those quantum field equations and using the fact that C, 
and rj have vanishing expectation values, we obtain a pair of equations equivalent to 








s + lQ{e+3v') + z{i^) - 



9 



4v^ 

9 
4^2 



(40) 
(41) 
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We will refer to these as the tadpole equations for ^ and for t], respectively. The quantum 
field equations for ^ and rj that correspond to ( [T^ ) are 



= e + (V^ + F)r/ + + ^Q[e + 3v-{e + 3v') 



+ V2Z - i^v)) - f {vie + V') - {vie + V'))) , (42) 



+ V2Q - m) - f {ae + v') - iae + v'))) ■ (43) 



The expressions (jS^) and (|21|) for the number density and the current density can be 
written as 

P = v^ + w' + ^{e + v'), (44) 

j = 2{vVw - wVv) + {^Vr]-r]VO- (45) 



Our strategy will be to solve the quantum field equations (^21) and (|^), treating the terms 
in ( |32D as perturbations, and then to calculate the expectation values in pO|), (PD, 
and ( ^51) to first order in ^fA. 



4 Semiclassical approximation 

For the perturbation theory defined by the decomposition of the action (|3l| ) into free and 
interaction parts, the appropriate dimensionless expansion parameter is the product of the 
coupling constant g and the coherence scale A. The semiclassical approximation is defined 
by truncating the expansion in powers of gA at first order. If the dimensionless ratio Aj ^fgp 
is held fixed, this is equivalent to truncating at first order in \fpc? . In terms of Feynman 
diagrams, accuracy to first order in gA corresponds to keeping only the contributions from 
one-loop diagrams. The expectation values in ( ^OD and (^Tj) that are cubic in ^ and t] receive 
contributions only from diagrams with 2 or more loops and therefore can be omitted. Thus 
the tadpole equations reduce to 

= T + iz(3e' + r/2) + Ql^in), (46) 

= ^ + ^Q(e' + 3r/2) ^ ^47) 

We will refer to these equations as the semiclassical tadpole equations for ^ and r/, respec- 



tively. In the quantum field equations (^21) and (^3]), the terms that are quadratic or cubic 
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in ^ and 7] contribute to the expectation values in ( ^61) and (^T]) only through diagrams 
with 2 or more loops. We therefore need to keep only the linear terms in the quantum field 
equations: 








-1) + (V' - A' + X)e + Rv- 



(48) 
(49) 



These equations describe the propagation of the fields ^ and r] in the presence of the sources 
X, Y, and R, but with no other interactions. Note that the explicit dependence of (|^) on 
A is cancelled by the A-dependence of the source X. We will find that a judicious choice 
of A will allow the source X to be treated as a perturbation, along with Y and R. The 
expectation values that appear in the semiclassical tadpole equations (^) and ( ^Tf ) and 
in the expressions for p and j in ( ^^ and ( ^5D are functionals of the sources X, Y, and 
R. Under the time-reversal symmetry (|2^), ^, X, and Y are even, while rj and R are odd. 



The expectation values (^^) and (r^^) are therefore even functionals of R, while (^77) and 
(^V?7 — ?7V^) are odd functionals of R. 



Since the sources X, Y, and R depend on v and w, the expectation values in (^4D-(^) 
are functionals of v and w. In order to obtain an equation for p and j, we will use (^4|) and 
(^5[) to eliminate v and w from (|46|) and (|47|) in favor of p and j. The expectation values 
that appear in these equations are the corrections to the classical equations from quantum 
field fiuctuations. They are smaller than individual terms in the classical equations by a 
factor of -s/po?. The semiclassical approximation includes all terms through first order in 
^/pa?. In the course of manipulating the equations (p4D-(^) to eliminate v and w, there is 
no loss of accuracy if we expand the resulting expressions to first order in the expectation 
values, dropping higher order terms. The resulting equations still includes all effects of 
quantum field fiuctuations through first order in y/pa-^. Solving (|4^) and (|45|) for v and 
Vw and expanding to first order in the expectation values, we obtain 

V = ip-wY^' - \.... {e + ^'), (50) 



Vw 



4(p-U72)l/2 



2v 



2v 

We can derive analogous expressions for Vf , V^f , and V^w by differentiating (^) and 
^TT). By substituting iteratively and expanding to first order in the expectation values. 



(51) 



we can eliminate v, Vf , V f , Vw, and V w from the right side of the equation in favor of 
p, j, their derivatives, and also w. The resulting expressions at a specific point fq become 
particularly simple if we use a global phase transformation to set w(ro) = 0: 



Vt;(ro) 



4pV2 
1 



ip' 



,3/2 



4p - {e + v') 



ApVp + Vp{e + v') - 2pV{e + v') 



(52) 
(53) 
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w{ro) 
Vw(ro) 



0, 



,3/2 



(54) 
(55) 

(56) 



16p5/2 

+ 4pVp- V(e2 + r72) - 4p2v'(e' + r/') + 8pj • (^Vr/ - r^V^ 



^P' 



3/2 



4pV ■ j + V ■ + f]') - 4pV ■ {^Vv - r?VO 



(57) 



Inserting the expressions (|52D-(|57D into ( ^61) and ( ^Tf ) and expanding to first order in 
the expectation values, we obtain expressions for the tadpole equations at the point fq 
with V and w ehminated in favor of p and j. The semiclassical tadpole equation ( ^7[ ) for rj 
reads 







8p3/2 



4pV ■ j + V • j + r/') - 4(7p2 (^r/) - 4pV ■ {^Vr] - r/VO 



(58) 



The last term in (|58|) can be simplified by using the quantum field equations (^8]) and 
for ^ and 77: 



(59) 



In the last step, we have used the expressions (pSf), (p6D, and (^) for the sources and the 
fact that + yf) is time- independent. Since w vanishes at the point yq and p — f ^ is of 
order -y/po^ at that point, (|59|) reduces to —gp{C,Ti), up to terms that are second order in 
y/pa^. Thus the last two terms in (|58|) cancel. Up to an overall multiplicative factor, the 
first two terms in (^) reduce to the time-independent continuity equation V ■ j = 0. 

The semiclassical tadpole equation (H) for at the point tq reduces to 



4p2 



^1/2 



4p5/2 



-pVV + (Vp)^+f (e^ + ry^) 



- ^ (p - - f p - ^ [-2pvv + ivpr+f]) {e + .^). 



(60) 
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The last line in (|60|) is proportional to the classical tadpole equation, which is just the 
first line of (pO). The classical equation differs from zero by terms of order \/pa^. Since 
the expectation value in the last line of (^) is also of order \/pa-^, the last line is second 
order in \/pa^ and can be omitted. Alternatively, it can be cancelled by multiplying (|60|) 
by 1 + (^^ + r/^) / (4p) and keeping only those terms that are first order in the expectation 
values. To complete the derivation of our equations for p and j, we must calculate the 
expectation values that appear in the remaining terms of (|60D and express them in terms 
of p and j. 



5 Thomas- Fermi expansion 

In the previous section, the semiclassical tadpole equation for ^ was evaluated at a point 
To where w vanishes and expressed in terms of p, j, and expectation values of operators. In 
the resulting equation (|60D, the operators constructed out of quantum fields C, and r] that 
satisfy (^) and (^D]). In this section, we use the Thomas- Fermi expansion to evaluate the 
expectation values and reduce (|60D to a partial differential equation involving p and j only. 



The expectation values in (^) are functionals of the sources X, Y, and R given in 



(|35D , (0), and (pQ]). We would like to expand the expectation values in powers of the 
sources and their gradients. The length scale for significant variations in the sources is 
much greater than the coherence length {gp)~^^'^- The expansion of the expectation values 
in powers of gradients of the sources is possible if the expectation values receive significant 
contributions only from modes of the quantum fields ^ and rj that have wavelengths of order 
((7p)~^/^ or less. This can be guaranteed by imposing an infrared cutoff that eliminates 
the contribution from modes with much longer wavelengths. A dependence on the infrared 
cutoff indicates a sensitivity to length scales much greater than the coherence length. We 
will find that the dependence on the infrared cutoff cancels when (|60|) is expressed in terms 
of p and j only. 

Our infrared cutoff allows an expectation value at the point r to be expanded in powers 
of gradients of the sources at the point r, with coefficients that are functions of X, Y, and 
R at the point r. A further expansion in powers of X, Y, and R is possible if these sources 
are at least first order in either the gradient expansion or in \fpa?. This is certainly not 
true in general, but we can make it true at a specific point fq by a judicious choice of 
the arbitrary parameter A. Note that the source R in ( |39D vanishes at the point Vq if the 
f/(l) symmetry has been used to set wijo) = 0. The sources X and Y in (|35|) and (|36|) do 



not appear to be higher order in the gradient expansion or in \fpa^. However, since the 
expectation values are already of order ^/pa^ and we are keeping only terms to first order 
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in \/pa^, we can use the classical equation T = to simplify these sources: 

X = A^-gv^--V%, (61) 

V 

Y = -gw^--Vh. (62) 

V 

From (|62D, we see that l^(ro) is second order in the gradient expansion, since w vanishes 
at the point Tq. We can arrange that X(ro) also be second order in the gradient expansion 
by a judicious choice of A. A convenient choice is 

A' = 9p{ro)- (63) 



Any other choice that differs from (|6^) by terms that are higher order in the gradient 
expansion or in y/ pa^ will lead to the same final equation for p and j. 

Having imposed an infrared cutoff and chosen a specific value for A, expectation values 
can be expanded in powers of X, Y, R, and their derivatives in a neighborhood of the 
point tq: 

(^^) = ao + aiX + a2Y + a^V^X + aeV^Y 

+ a^X^ + a^iVXf + ajR^ + asiVR)'^ + (64) 
= bo + hX + hY + fogV'X + beV^Y 

+ b^X^ + b^iyXf + hR^ + bs{VRf + (65) 
(e^) = ciV'i? + CaVX-Vi? + (66) 
(^Vr^-r/VO = diWR + .... (67) 

The coefficients a^, 6j, q, and di are functions of A. They can be obtained by calculat- 
ing Feynman diagrams using the methods illustrated in Appendix p. These calculations 
are streamlined by using dimensional regularization to control the infrared and ultraviolet 
divergences that appear in the individual diagrams. Analytic expressions for the dimen- 
sionally regularized integrals are given in Appendix 0. The final results for the coefficients 
are given in Appendix |y. We have explicitly shown only those terms in (p^)-(p7|) that are 
required to calculate the condensate equation, the number density, and the current density 
through second order in the gradient expansion. Because of time reversal symmetry, (^^) 
and {7f) contain only terms that are even in R, while (^77) and (^V?7 — rjVC) contain only 
terms that are odd in R. Differentiating the sum of (|6^ and (^) and keeping only those 
terms that contribute through second order in the gradient expansion at the point fq, we 
obtain 

V(e' + r/') = (ai + 6i)VX, (68) 
V'(e' + ?]') = (ai + 6i)V'X + (a2 + 62)W 

+ 2(a4 + 64)(VX)=^ + 2{a7 + h){\7Rf. (69) 



17 



We proceed to express X, Y, and R and their derivatives at the point Tq in terms 
of p and j. After differentiating (|HT|), (p^), and ( |50D and evaluating them at tq, we can 
use (^)-([57|) to ehminate v, w, and their derivatives. Since X, Y, and R appear only in 



expectation values that are already of order ypcfi, the equations (M)~(|57|) can be simplified 
by dropping the expectation values: 



Vt;(ro) 
V«;(ro) 



P 



1/2 



4p3/2 



-2pVV+(Vp)2+f 



0, 



2pi/2 

1 



2pi/ 



The resulting expressions for the sources X, Y and and their derivatives are 

VX(ro 
V2X(ro 



vr(ro 

W(ro 

i?(ro 
Vi?(ro 

V'i?(ro 



(70) 
(71) 

(72) 
(73) 
(74) 

(75) 



1 

V 


-2pVV+(Vp)2+f 


(76) 


= -^?vp, 




(77) 


^ ■ 

2p^ 


-2pVV+j' 




(78) 


1 

V 


-2pVV+(Vp)2+f 


(79) 


= 0, 






(80) 


9 
2p 






(81) 


= 0, 






(82) 


^ . 




(83) 


2p 


[Vp-j + pV 


■j]- 


(84) 



Inserting the expressions (|6^ -(|69D for the expectation values into (^), using the expres- 
sions ([r^)-(|5^) for the sources and their derivatives, and truncating at second order in the 
gradient expansion, the semiclassical tadpole equation for ^ at the point rg reduces to 







4p2 



-2pVV + (Vp)2+f 



.1/2 



(ao + ho) + (2ai + 02 + &i)A^ + 2a3A^ 



4p3/2 
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2(ao + bo) + (3ai + 03 + 2bi)A^ + 4(04 + &4)A^ + 4a5A^ 
2(ao + 60) + (2ai + 61 - 62 + Sd^A^ 
+ (2a3 - 2a6 + 07 + 67) A^ + agA^ 



,5/2 



We have omitted the last term in (^), because it is second order in y^po^. 



^5) 



Inserting the coefficients Oj, hi, and di from Appendix |C| into (^), using the expression 
for A^ in and multiplying by — p^^/^, the semclassical tadpole equation reduces to 



P 



V + y + 



4p2 



+ 



127r- 



-2pVV+(Vp)2+f 
17 



192p2 



4pV> + (Vp)^ 



^6) 



This is an algebraic relation between the values of p, j, and derivatives of p. It was 
derived at a specific point ro where the condensate is real-valued. However a U{1) 
symmetry transformation can be used to make real-valued at any given point r, and p 
and j are independent of that phase. Therefore (^) must be valid at any point where 
the semiclassical and Thomas-Fermi approximations can be justified. Using dimensional 
analysis to reinstate the appropriate factors of h and 2m and using (|^) to express g in 
terms of the scattering length, we obtain our final result (|^). 

The coefficients Qq and 60 in (HD are cubicly ultraviolet divergent, while Oi, 02, 61, &2, 
and di are hnearly ultraviolet divergent. These power ultraviolet divergences are set to 
zero by dimensional regularization. If we had used a momentum cutoff A(/v, the ultraviolet 
divergences would have to be cancelled explicitly by counterterms. If those counterterms 
were chosen to be purely cubic or linear in Auv, then the coefficients after renormalization 
would be identical to those obtained directly using dimensional regularization. Thus dimen- 
sional regularization systematicly throws away power divergences that would ultimately be 
cancelled by counterterms. 

The coefficients a^, ag, 62 and 67 in the quantum correction proportional to are 
logarithmicly infrared divergent, and must be calculated using an infrared cutoffi With 
dimensional regularization, these logarithmic divergences appear as poles in D~3, where D 
is the number of spatial dimensions. However these divergences cancel in the combinations 
of coefficients that appear in (^51). The cancellation of the infrared divergences indicates 
that our final equation ( p6D is insensitive to length scales much greater than the coherence 
length. To obtain the cancellation of infrared divergences, it was crucial to have omitted 
the last term in (pDj). That term is infrared divergent, but the classical equations can be 

2^. The infrared divergence can be cancelled only 



used to show that it is second order in 
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by other terms that are second order in \/pa^. All such terms have been systematically 
dropped. 

There is a remarkable cancellation of the correction proportional to {g^ pY^'^y^ / in 
(|85D . Because of this cancellation, the current j does not appear at all in the semiclassical 
correction term in (|86|). We have no explanation for this cancellation. 



The semiclassical correction term in (^) is a functional of p and j that includes a 
multiplicative factor of {g^pY^'^. The Thomas- Fermi approximation has been used to ex- 
pand the cofactor of {g^pY^"^ to second order in ^/-R, where ^ = (gp)'^^"^ is the coherence 
length and R is the length scale for significant variations in p. Since ^ scales like p~^^'^, the 
Thomas-Fermi expansion breaks down in regions where p approaches zero, such as near 
the edge of the condensate or near the core of a vortex. In such a region, the cofactor of 
(g^pY^^ in ([86| ) becomes a nonlocal functional that depends on the values of p and j at 



points within a distance of order R. However the multiplicative factor of {g^pY^"^ makes the 
semiclassical correction negligible in regions where p — >■ 0. If ((7^p)^/^/(127r^) is sufficiently 
small, then each term in the semiclassical correction in (|86[) will be small compared to one 
of the terms in the classical equation. It therefore does no harm to include these terms. 
Thus the equation ( ^61) can be used everywhere, even near the edge of the condensate or 
near the core of a vortex. 

The chemical potential p is the change in the total energy if a single atom is added 
to the condensate. The terms on the right side of (|86D therefore have simple physical 
interpretations as contributions to the energy of that additional atom. In the mean-field 
approximation, that energy consists of the potential energy V , the interaction energy gp/2, 
and the gradient energy. The signs of the semiclassical corrections are such as to increase 
the interaction energy and decrease the gradient energy. 



6 Condensate equation 

In this section, we calculate the expectation values in the semiclassical equation ([15|) for the 
condensate to second order in the gradients of using an infrared cutoff. The condensate 
equation is expressed in the form of a partial differential equation for that depends 
logarithmicly on the infrared cutoff. We will find that if is eliminated in favor of p and 
j, the dependence on the infrared cutoff cancels and we recover the equation (|86D. 



The expectation values appearing in the condensate equation ( ]T3| ) are expectation 



values of operators involving real quantum fields and 77 that satisfy (H) and (|49D: 

{^H) = lie + v'), (87) 
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(88) 



Expansions for (^^), {rf), and (^77) in powers of X, F, i?, and their derivatives are given in 
(|6^-(|66D. These expansions are vahd in a neighborhood of a point fq where w vanishes, 
provided that we make a suitable choice for the arbitrary parameter A. If we wish to express 
the expectation values (|87|) and (|^) in terms of the condensate, the simplest choice is 



^9) 



This differs from the choice ( |63D made in the previous section by terms of order \J pa-' 
The choice (^Qf ) leads to the same final result if we truncate at the same orders in \/pc? 
and in ^/R. After differentiating the expressions (pTj), (p^), and ( P^ for the sources and 
evaluating them at the point tq, they can be simplified using w(ro) = 0. Inserting the 
resulting expressions into (|64|)-(|6^), we obtain 



ao - 

- 2 
bo - 

- 2 



ai + 02 + 2a3A^ V^v/v 
61 + 62 + 263A2I V\/v 



CigvV'^w — 2 
digvVw. 



g{Vv) 



2ae-asA^ g{Vw)\ 
giywf, 



2h 



ci — 02^^ gVv ■ Vw, 



(90) 

(91) 
(92) 
(93) 



Some of the coefficients in the expansions (^O]) and (|9l|) are infrared divergent and must 
be calculated using an infrared cutoff. The coefficients and ag are logarithmicly infrared 
divergent, but the divergence cancels in the combination 2a6 — A^ag that appears in (|90|). 
The coefficients and 6g are quadraticly infrared divergent, but the leading divergences 
cancel in the combination 266 — A^6g that appears in (|9T[), leaving a logarithmic divergence. 



There are also logarithmic infrared divergences in 62, 63, 65, and 67. We choose to use 
dimensional regularization for the infrared cutoff. Logarithmic infrared divergences appear 
as poles in Z) — 3, where D is the number of dimensions. The expressions for the infrared 
divergent integrals can be simplified by trading the poles in D — 3 for logarithms of a 
momentum scale piR defined by 



1 

log — 



D-3 



, 4A2 

+ log — 2 + ^' 



(94) 



where 7 is Euler's constant and /i is the momentum scale introduced by dimensional regular- 
ization. The momentum scale fim can be interpreted as a conventional infrared momentum 
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cutoff. Using tlie values for the coefficients given in Appendix |C| and using tlie cfioice (|8 
for A, tlie expectation values (PD|)-(P^ reduce to 



2 



144vr2 
1 



(95) 



127r 
5 



967r2 



log — 



647r2 
2 44 



log — 



62 
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^1/2 (v^t; + (Vu;) V 



g'/\Vvr/v 



1 

— (V^U7 + ■ Vw/v 



127r2 



127r 



(96) 

(97) 
(98) 



The infrared-divergent coefficients in ( P^ ) indicate that (?7^) is logarithmicly sensitive to 
length scales much greater than the coherence length. 



Inserting (|95|)-(p7D into (|87| ) and (|88D , we obtain expressions for the noncondensate 
density and the anomalous density that hold at a particular point Tq where w{ro) vanishes. 
We can use the U{1) symmetry to deduce general expressions for these expectation values 
in terms of and gradients of 0. The noncondensate density ( |57D is invariant under 
the U{1) symmetry. At leading order in the gradient expansion, it must be a function 
of |0p = 0*0, which reduces to at the point tq. Thus the f/(l)-invariant expression 



for the V term in ([87| ) is simply 



At second order in the gradient expansion, the 



noncondensate density can be expressed as a linear combination of V0* ■ V0, V^(0*0), 
and [V(0*0)]^, with coefficients that are functions of |0|. Any f/(l)-invariant expression 
that is second order in gradients of can be expressed as a linear combination of these 
three terms. We can deduce the coefficients of those terms in (ip^ip) by evaluating the 
linear combination at the point Fq and comparing to (|87D. The resulting [/(l)-invariant 
expression for the noncondensate density is 



1 



247r2- 



+ 



^3/2 I 



7 



3847r2 



log 



15367r2 



log — 7^ 4 



g 



1/2 1 



22 

[V(0* 



iV2| 



-1^2 



(99) 



The anomalous density (|88D is a functional of that transforms like ■0^ under the U{1) 
symmetry (^Bj). At leading order in the gradient expansion, it must have the form 0^ 
multiplied by a function of |0|. Thus the f/(l)-covariant expression for the term in ( |8^ ) 



IS 



At second order in the gradient expansion, there are five independent terms that 
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transform like ijj'^ under the f/(l) symmetry: 0V^0, (V0)^, 0^V0* ■ V0, 0^V^(0*0), and 
0^ [V(0*0)]^. The anomalous density ( |5^ ) must be a linear combination of these terms with 
coefficients that are functions of |0|. We can determine the coefficients by evaluating this 



linear combination at the point Tq and matching with 
expression for the anomalous density is 



I). The resulting f/(l)-covariant 



87r2 



3/5 



3847r2 
7 

1536^ 



log 



247r2 

26' 



log 



15 
64 
" 21 



g'/'\<P\-V [V(0^ 



(100) 



Inserting the expectation values ( pQ]) and (|100|) into (p!5D, we obtain the semiclassical 
equation for the condensate to second order in the Thomas-Fermi expansion: 







1 



1 



5 



^9 



3/2 I 



+ 



+ 



+ 



487r2 
5 

7687r2 
7 

30727r2 



247r2 

U*(V0)2 + 30V0* ■ V0 



1 + 



247r2 



(7|0|2 158 
log ^ TF 

^ ^ f^h. 21 , 



(101) 



The infrared-divergent coefficients in the last two terms of ( p.01| ) indicate that, at second 
order in the gradient expansion, the condensate is logarithmicly sensitive to length scales 
much greater than the coherence length. 

We can also obtain f/(l)-invariant expressions for the number density and the current 
density in terms of the condensate (p. For the number density, this is simply a matter of 
inserting the noncondensate density (p^) into (EOf): 



+ 



247r2 
7 



15367r2 



1 9 

log -2 



38471^ 



, 9\<P\' 
log^— 



22 



,V2| 



-4 



(102) 



To obtain the current density, we need the [/(l)-invariant expression for the expectation 
value (|9|): 



3/2 1 



247r2^ 
23 



*V0 - 0V0*) 



(103) 



Inserting this into (plD, we find that the current reduces to 

1 



—I 



247r2 



9"'\<P\ 



- 



(104) 



Given the condensate equation ( |101| ) and the expressions for p and j in (|102|) and 
( |104| ), it is straightforward to derive the partial differential equations for p and j given in 

and (0). It is convenient to first multiply ( |101|) by 0*. The imaginary part of the 
resulting equation immediately gives the continuity equation V ■ j = 0, where j is given 
in (|104|) . In the real part of the equation, we can use (|102|) and (|104|) to eliminate in 
favor of p and j. We find that all the infrared divergences can be organized into a term 
that is proportional to the classical equation for p and j. That term is second order in 
\fp(? and can be consistently deleted. The remaining terms reduce to (|86D up to an overall 
multiplicative factor. The fact that the [/(l)-invariant form of the condensate equation 
(|101|) reproduces the differential equation (|86D for p and j provides a stringent check on 
both results. 



7 Conclusions 

In this paper, we have calculated the corrections from quantum field fluctuations to the 
mean-field equations for a time-independent state of a Bose-Einstein condensate. We used 
only controlled approximations that correspond to truncations of systematic expansions in 
small quantities: the semiclassical approximation, which is the expansion to first order in 
y/pa"^, and the Thomas- Fermi approximation, which is an expansion in powers of ^/R. The 
consistency of these approximations is guaranteed by the fact that they are controlled by 
small expansion parameters. 

By integrating out the quantum field fluctuations, we obtained remarkably simple self- 
contained equations for the condensate (p and for the hydrodynamic fields p and j. The 
effects of quantum field fluctuations were taken into account through local correction terms 
to the partial differential equations of the mean-field approximation. In the Thomas-Fermi 
limit, the correction to the Gross-Pitaevskii equation is the |0|^0 term in (^). At second 
order in the Thomas-Fermi expansion, the corrections to the Gross-Pitaevskii equation can 
no longer be expressed in a local form. However the effects of quantum field fluctuations 
can be taken into account by adding local correction terms to the hydrodynamic form of 
the mean-field equations. The resulting equations are the continuity equation (|^) and the 
nontrivial equation (^). The numerical solution of these partial differential equations is 
much simpler than solving the equations of the Hartree-Fock-Bogoliubov approximation. 

The Thomas-Fermi approximation has also been used in a recently developed varia- 
tional Thomas- Fermi theory of a nonuniform Bose-Einstein condensate |]16|. In the case 
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of atoms that interact only through the S-wave scattering length, this approach leads to 
results similar to ours in the Thomas- Fermi limit. The authors were unable to proceed 
beyond leading order in the Thomas-Fermi expansion because of the breakdown of the 
gradient expansion. We avoided this problem by changing to the hydrodynamic variables 
p and j. 

Our methods can be readily extended to time- dependent states of a Bose-Einstein 
condensate. The resulting equations would describe the collective excitations of a Bose- 
Einstein condensate, including the leading effects of quantum field fluctuations. A partic- 
ularly simple application of these equations would be to study the frequencies for small 
oscillations of the condensate. In the mean-field approximation and in the Thomas-Fermi 
limit, the oscillation frequencies for a spherically-symmetric harmonic oscillator potential 
V{r) = are particularly simple The frequencies are independent of the scat- 



tering length a and ratios of frequencies are independent of cu. Quantum field fluctuations 
give fractional corrections to the frequencies proportional to N^/^{a?muj /fiY^^ . It would 
be very interesting to calculate these changes in the frequencies and see if they can be 
measured in experiments. 

Up to this point, we have only applied our methods to the relatively simple problem 
of a Bose-Einstein condensate at zero temperature. It may be possible to generalize our 
methods to the case of nonzero temperatures T that are below the critical temperature 
Tc for the Bose-Einstein phase transition. The temperature introduces a new length scale 
into the problem, the thermal wavelength h / y/2'mkT, and therefore a new dimensionless 
ratio mkT / (h^ pqo) . The condensate number density po decreases to as T approaches 
Tf.. Since the coherence length ^ scales like (poc^)"^^^, the Thomas- Fermi expansion in 
^/R breaks down for large enough T. However the product of \/pqc? and mkT / {h^ p^a) 
remains small at all temperatures below the critical region near Tc [|1^]. This quantity can 



therefore be used as an expansion parameter for a controlled approximation. One may 
be able to use this controlled approximation to systematically integrate out the effects of 
short-wavelength fluctuations around the mean field. It would be interesting to see if this 
can lead to a simplification of the equations for nonzero temperature that is comparable 
to the simplification that we have achieved at zero temperature. 
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A Integrals 



In this Appendix, we give analytic expressions for the energy integrals and the momentum 
integrals that are required to calculate the coefficients in the expansions of the expecta- 
tion values (|6^-(|67D in powers of the sources and their derivatives. We use dimensional 
regularization to cut off any infrared or ultraviolet divergences in the momentum integrals. 

The energy integrals have poles at u; = ±e, where e = k^/h?~+~K^. The poles come 
from Feynman propagators, which are defined by an iO'^ prescription. The integrals can 
be evaluated using contour integration. The specific integrals that are required are 



where {z)n is the Pochhammer symbol: {z)n = ^{z + n)/r{z). 

Some of the momentum integrals are ultraviolet divergent or infrared divergent or both. 
We choose to use dimensional regularization to regularize both ultraviolet and infrared 
divergences. This involves calculating the integral as an analytic function of the number 
of dimensions D and analytically continuing it to Z) = 3. The analytic continuation sets 
power ultraviolet divergences and power infrared divergences equal to 0, but logarithmic 
divergences appear as poles in D — 3. It is convenient to take the integration measure to 
be fi^~^ J d^k/{27c)^ , where fi is an arbitrary momentum scale. The prefactor fi^~^ has 
been inserted so that the regularized integral has the same engineering dimensions as in 



The insertion of a source into a loop diagram shifts the loop momentum k by the 
momentum flowing through the source. The gradient expansion corresponds to expanding 
the loop integral in powers of the momenta flowing through the sources. This expansion 
generates integrals over the momentum k with vector indices. These integrals can be 
reduced to scalar integrals by averaging over angles in D dimensions using the formulas 




(Al) 



(A2) 



D 



3. 




(A3) 



(A4) 



One must avoid setting D = 3 on the right side of 
logarithmic divergence. 



( |A3| ) if the scalar integral over k has a 
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The scalar momentum integrals can be written in the form 



/i 



3-D 



2\m 



(A5) 



where m and n are integers. Evaluating ( |A5|) in a region of D where it is convergent and 
then expressing it as an analytical function of D, we obtain 



'■m,n 



T{n-m- f )r(i^ + m) , 3.^ 



(4vr)«/2r(f)r,, 



These integrals satisfy the identities 



(A6) 



T — T — r 

-'m— l,n— 2 -'m+l,n ■'^ ^m,ni 

[D + 2m- n)Im,n = nIm+2,n+2- 



(A7) 
(A8) 



The first is just an algebraic relation, while the second follows from applying integration 



by parts to ( |A5| ) 



One of the advantages of dimensional regularization is that the formulas ( [A^ 
hold even if the integrals are infrared or ultraviolet divergent. For D = 3, Im,n is ultravio- 
let divergent if m and n satisfy m — n > — |, and infrared divergent if 2m — n < —3. The 
ultraviolet divergences are power divergences that are set to zero by dimensional regular- 
ization. Thus, unless 2m — n < —3, we can set D = 3 in ([A6|) . The integral then reduces 
to 



r(n-m-|)r(^ + m) 
47r2r(f) 



n < 2m + 3 . 



(A9) 



If 2m — n = —3, Im,n has a logarithmic infrared divergent for D = 3. With dimensional 
regularization, this divergence appears as a pole in D — 3. After extracting the pole from 
the factor r(-^^) in ( |X6| ) and expanding the remaining factors to first order in D — 3, we 
obtain 



1-2,-1 
^-1,1 
-^0,3 



1 



47r2 
1 

4^ 
1 

4^ 



2 , 4A2 

+ log — + 7 - 4 



D-3 
2 

D-3 
2 

D-3 



, 4A2 

log— -+7-2 



+ log 



4A2 



+ 7-4 



A, 
1 

A' 
1 

A3' 



(AlO) 
(All) 
(A12) 



If 2m — n = —5, Im,n has a quadratic infrared divergence for D = 3, underneath which is a 
logarithmic infrared divergence. Dimensional regularization sets the quadratic divergence 
to zero, but the logarithmic divergence appears as a pole in D — 3. After extracting the 
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_D-5^ 



pole from the factor ^ 
D — 3, we obtain 



-3,-1 



-2,1 



m 



and expanding the remaining factors to first order in 



D-3 



4A2 
log — ^ + 7 



87r2 



4A2 



1 

A' 



D-3 



log- 



+ 7-3 



1 

A3- 



(A13) 
(A14) 



The product of -D — 3 with a pole in D — 3 gives a nonzero result in the limit D — > 3. Thus 
if any of the infrared divergent integrals in (|A10|) -( [A14| ) is multiplied by a function of D, 
one must take care to expand that function to first order in D — 3 before taking the limit 
D^3. 



B Feynman diagrams 

In this Appendix, we illustrate the use of Feynman diagrams to calculate the coefficients 
in (|6^-(|67|). In these equations, the expectation values of quantum field operators in the 
presence of the sources X, Y, and R are expanded in powers of the sources and their 
gradients. As an example, we will calculate the coefficient Ci in the expansion (^) for 

The basic ingredients of Feynman diagrams are propagators and vertices. It is con- 
venient to calculate the diagrams in momentum space. The diagonal propagators for the 
fields ^ and t] are represented by solid lines and dashed lines, respectively. Their Feynman 
rules are 

(Bl) 



w2-e2(/t)+z0+' 



(B2) 



where k is the momentum and uj is the energy fiowing through the line. The off-diagonal 
propagator for ^ and rj is represented by a line that is half solid and half dashed. Its 
Feynman rule is 

c^2-e2(^)+^0+' ^^^^ 

where uj is the energy flowing from the dashed end of the line toward the solid end. Each 
of the operators i^^, rf, and ^i], and ^Vrj — rj'V^ creates two outgoing lines. The vertex 
factors for the operators 77^, and ^77 are 2, 2, and 1, respectively. The vertex factor for 
the operator ^Vr] — f]VC, is p — 2k if the momenta flowing out the ^ and r] lines created 
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by the operator are k and p — k, respectively. The only other vertices needed to calculate 
the one-loop diagrams for the expectation values of these operators are for the external 
sources X, Y, and R. If momentum p is flowing out of the diagram though one of those 
sources, the vertex factor is iX{p), iY{p), or iR{p). The one-loop diagrams are integrated 
over the energy u and momentum k running around the loop, with the measures / duj/{2'K) 
and / d^k/{2'n')^, where D is the number of spatial dimensions. Finally, there is a 
symmetry factor of 1/2 if the diagram has a reflection symmetry. 

We consider the diagrams in Fig. |l|, which represent the contributions to the expectation 
value {^rj) involving one insertion of the source R. If the external momentum p enters 
through the operator ^r/ and exits through the source -R(p), the expression for the first 
diagram in Fig. ffl is 



2^^ 



3-D 



-UJ 



(27r)^C^2_,2(|k + p|) 



zR{p)- 



e^(|k|)- 



(B4) 



We have written the Feynman rules for each of the propagators and the vertex in the loop 
in the order in which they appear as you go counterclockwise around the loop. There is 
an implied +zO"'" prescription in the denominator of each of the propagators. 

The first step in evaluating the diagram is to expand the integrand in powers of the 
external momentum p. The expansion of the first denominator in (p^) to second order in 
p has the form 



e^(|k + p|) 



1 



+ 



(^' + ^) (p' + 2p-k) + 4(p-k)^ 



+ 4 + 



fc2 



(P ■ k)^ 



1 



(Cj2 - e2)3 



/1\2 



(B5) 



where e = e(fc) on the right hand side. We can average over the angles of k by substituting 
kikj —>■ k'^6ij/D and fcj 0. The term linear in p ■ k drops out, and the diagram ( P^ ) 
reduces to 



iR{p) 



duo 
2^ 



3-D 



d^k 



(27r)« (cu2 



c2^2 



1 



We can now use the formula (^) to evaluate the integrals over uj. This reduces the diagram 
to an integral over k: 



R{p)^ 



3-D 



(27r)^ \Ae 



+ 



D{k^ + e')+Ak^ {k^ + e''f 
IQDk^e^ ^ ~8Dkh^ 



P 



(B7) 
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In terms of the integrals Im,n defined in 



this reads 



r^o,i^(p) + 



1 1 D-2 

1 Q c; — 1 1 Q — 

8D ' 16 ' IQD 



/-I, 



(B8) 



The expression for the second diagram in Fig. |l] is 

duj 



27T 



3-D 



It is calculated in the same manner and the result is 



r4,i^(p) + 



8D 



'3,5 



"-'1,3 + — TTTT^-'-l,! 



16D 



p2i?(p). 



The sum of the diagrams ( [BSD ^'^d ( |B1U| ) is 



1 



h,iR{v) + 



'3,5 



1 

AD 



-h 



,3 



p2i?(p) 



Using the identity 



this reduces to 
1 



1 



/0,li?(p) - ^/l,3p'i?(p). 



After Fourier transforming to coordinate space, we can write this as 
(er?(x)) 



(B9) 



(BIO) 



(Bll) 



(B12) 



(B13) 



where . . . represents term with more gradients of R or more powers of X, F, or R. Com- 
paring with the expansion for (^77) in (|66|), we can read off the coefficient ci = /i 3/(2D). 
We also find that there is a term cqR in the expansion (pBD, with cq = -^o,i/2. This term 
was omitted in (|66D, because cq does not contribute to any of the quantities calculated in 
this paper. 



C Coefficients 



In this Appendix, we express the coefficients that appear in (p^)-(p7D in terms of scalar 
momentum integrals. The coefficients can be calculated by evaluating Feynman diagrams, 
as described in Appendix ^ of Ref. |J7[. After using ( |^ ) and ( [A2D to integrate over the 
energy and ( |A3| ) and ( |A4| ) to average over angles, the momentum integrals can be reduced 
to the scalar integrals Im,n defined in Appendix 0. 
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There are twelve coefficients that are required to calculate and p for the ground state 
in the semiclassical approximation and to second order in the Thomas-Fermi expansion. 
They were calculated previously in Ref. |0 using a large momentum cutoff to regularize 
ultraviolet divergences and a small momentum cutoff to regularize infrared divergences. 
Since we are using dimensional regularization for both infrared and ultraviolet divergences, 
we need the expressions for those coefficients in D dimensions. Eight of the coefficients are 
given by the same expressions as in Ref. [0: 



n ^ — 
Uq — 


2^1,1, 


(CA) 


CLi = 


4^2,3, 




a2 = 




(C3) 


a4 = 


Y^^3,5, 


(C4) 


bo = 


2^-1--!' 


(C5) 


bi = 




(C6) 


b2 = 


4^-2,-1, 


(C7) 


64 = 


16 • 


(C8) 



The expressions for the remaining four coefficients depend explicitly on D: 

as = Y^[5/5j- 2/3,5 + /i,3], (C9) 

«5 = 7^[35/6.9- 10/4,7 + 3/2,5], (CIO) 

h = -Y^[/3,5 + 2/1,3 + /-i,i], (Cll) 

= -^[^4,7 + 2/2,5 + /o,3]. (C12) 

The factor of D in the denominator arises from the angular average in ( |A3|) . We have 
used the identity (|A8|) to put the expressions (|C9|) -( |CT2D into a standard form with no 



factors of D in the numerator. Setting D = 3, we recover the expressions in Appendix C of 
Ref The integrals /-i,i and /o,3 are logarithmicly infrared divergent. With dimensional 
regularization, they have poles in D — 3. There are therefore nonzero contributions to 63 
and 65 that arise from expanding the prefactor 1/D in (|C11|) and (|C12| ) to ffist order in 
D-3. 
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There are nine additional coefficients tfiat are required to calculate 0, p, and j for a 
general time-independent state. The expressions for these coefficients in D dimensions are 




i-[/3,5 + 2/1,3 + 



(C13) 



(C14) 



(C15) 



(C16) 



(C17) 



(C18) 



(C19) 



(C20) 



^1 — 77h[-^2,3 - /o,l], 



(C21) 



The integrals /-3,-i and I-2,i are quadraticly infrared divergent, while and /o,3 are 
logarithmicly infrared divergent. With dimensional regularization, they have poles in /P — 3. 
There are therefore nonzero contributions to Oe, Og) ^6; and that arise from expanding 
the pref actor 1/D to first order in /P — 3. 
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FIGURE CAPTIONS: 



Figure 1: Feynman diagrams for calculating the coefficient ci.The energy and momentum 
flow counter-clockwise around the loop. 




Figure 1: Feynman diagrams for calculating the coefficient Ci. The energy and momentum 
flow counter-clockwise around the loop. 



34 



